﻿using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;

namespace SmartMathLibrary.DataMining.ClusterAnalysis
{
    /// <summary>
    /// This class provides a cluster analysis by using the CMeans algorithm. The implemented algrithm base on 
    /// the implementation from http://datamining.codeplex.com/.
    /// </summary>
    [Serializable]
    public class CMeansClusterAnalyser : AbstractClusterAnalyser
    {
        /// <summary>
        /// This field holds the clusters for the cluster analysis. 
        /// </summary>
        private readonly Cluster[] clusters;

        /// <summary>
        /// The needed iterations until the algorithm convergence.
        /// </summary>
        private int neededIterations;

        /// <summary>
        /// The fuzzyness factor to use.
        /// </summary>
        private double fuzzynessFactor;

        /// <summary>
        /// This field holds the data for the cluster analysis. 
        /// </summary>
        private ClusterAnalysisItem[] analysisData;

        /// <summary>
        /// This field holds the distance measure mode.
        /// </summary>
        private DistanceCalculationMode distanceMeasureMode;

        /// <summary>
        /// The objective function for the algorithm.
        /// </summary>
        private double objectiveFunction;

        /// <summary>
        /// This field holds the membership matrix.
        /// </summary>
        public double[,] membershipMatrix;

        /// <summary>
        /// The internal epsilon, which will be used by the algorithm.
        /// </summary>
        public const double Epsilon = 1e-5;

        /// <summary>
        /// Initializes a new instance of the <see cref="CMeansClusterAnalyser"/> class.
        /// </summary>
        /// <param name="analysisData">The data for the cluster analysis.</param>
        /// <param name="clusters">The clusters for the cluster analysis.</param>
        public CMeansClusterAnalyser(ClusterAnalysisItem[] analysisData, Cluster[] clusters)
            : this(analysisData, clusters, 2.0, DistanceCalculationMode.Euclidean)
        {
        }

        /// <summary>
        /// Initializes a new instance of the <see cref="CMeansClusterAnalyser"/> class.
        /// </summary>
        /// <param name="analysisData">The data for the cluster analysis.</param>
        /// <param name="clusters">The clusters for the cluster analysis.</param>
        /// <param name="fuzzynessFactor">The fuzzyness factor.</param>
        public CMeansClusterAnalyser(ClusterAnalysisItem[] analysisData, Cluster[] clusters, double fuzzynessFactor)
            : this(analysisData, clusters, fuzzynessFactor, DistanceCalculationMode.Euclidean)
        {
        }

        /// <summary>
        /// Initializes a new instance of the <see cref="CMeansClusterAnalyser"/> class.
        /// </summary>
        /// <param name="analysisData">The data for the cluster analysis.</param>
        /// <param name="clusters">The clusters for the cluster analysis.</param>
        /// <param name="fuzzynessFactor">The fuzzyness factor.</param>
        /// <param name="distanceMode">The distance measure mode to use during the cluster analysis.</param>
        public CMeansClusterAnalyser(ClusterAnalysisItem[] analysisData, Cluster[] clusters, double fuzzynessFactor,
                                     DistanceCalculationMode distanceMode)
        {
            if (analysisData == (ClusterAnalysisItem[]) null)
            {
                throw new ArgumentNullException("analysisData");
            }

            if (analysisData[0] == (object) null)
            {
                throw new ArgumentNullException("analysisData");
            }

            this.clusters = clusters;
            this.analysisData = analysisData;
            this.fuzzynessFactor = fuzzynessFactor;
            this.distanceMeasureMode = distanceMode;
            this.membershipMatrix = new double[this.analysisData.Length,this.clusters.Length];

            for (int i = 0; i < this.analysisData.Length; i++)
            {
                double sum = 0.0;
                ClusterAnalysisItem p = this.analysisData[i];

                for (int j = 0; j < this.clusters.Length; j++)
                {
                    Cluster c = this.clusters[j];

                    double diff = this.CalculateDistance(p, c);

                    if (diff == 0)
                    {
                        membershipMatrix[i, j] = Epsilon;
                    }
                    else
                    {
                        membershipMatrix[i, j] = diff;
                    }

                    sum += membershipMatrix[i, j];
                }

                double sum2 = 0.0;

                for (int j = 0; j < this.Clusters.Length; j++)
                {
                    membershipMatrix[i, j] = 1.0 / Math.Pow(membershipMatrix[i, j] / sum, 2.0 / (fuzzynessFactor - 1.0));
                    sum2 += membershipMatrix[i, j];
                }

                for (int j = 0; j < this.Clusters.Length; j++)
                {
                    membershipMatrix[i, j] = membershipMatrix[i, j] / sum2;
                }
            }

            this.RecalculateClusterIndexes();
        }

        /// <summary>
        /// Gets the clusters for the cluster analysis.
        /// </summary>
        /// <value>The clusters for the cluster analysis.</value>
        public Cluster[] Clusters
        {
            get { return clusters; }
        }

        /// <summary>
        /// Gets or sets the data for the cluster analysis.
        /// </summary>
        /// <value>The data for the cluster analysis.</value>
        public ClusterAnalysisItem[] AnalysisData
        {
            get { return analysisData; }
            set { analysisData = value; }
        }

        /// <summary>
        /// Gets or sets the distance measure mode to use.
        /// </summary>
        /// <value>The distance measure mode to use.</value>
        public DistanceCalculationMode DistanceMeasureMode
        {
            get { return distanceMeasureMode; }
            set { distanceMeasureMode = value; }
        }

        /// <summary>
        /// Gets the needed iterations of the algorithm.
        /// </summary>
        /// <value>The needed iterations of the algorithm.</value>
        public int NeededIterations
        {
            get { return neededIterations; }
        }

        /// <summary>
        /// Gets or sets the fuzzyness factor.
        /// </summary>
        /// <value>The fuzzyness factor.</value>
        public double FuzzynessFactor
        {
            get { return fuzzynessFactor; }
            set { fuzzynessFactor = value; }
        }

        /// <summary>
        /// Recalculates cluster indexes.
        /// </summary>
        private void RecalculateClusterIndexes()
        {
            for (int i = 0; i < this.analysisData.Length; i++)
            {
                double max = -1.0;
                ClusterAnalysisItem p = this.analysisData[i];

                for (int j = 0; j < this.clusters.Length; j++)
                {
                    if (max < membershipMatrix[i, j])
                    {
                        max = membershipMatrix[i, j];

                        if (max == 0.5)
                        {
                            p.ClusterIndex = 0.5;
                        }
                        else
                        {
                            p.ClusterIndex = j;
                        }
                    }
                }
            }
        }

        /// <summary>
        /// Perform one step of the algorithm.
        /// </summary>
        public void Step()
        {
            for (int c = 0; c < this.clusters.Length; c++)
            {
                for (int h = 0; h < this.analysisData.Length; h++)
                {
                    double top = CalculateDistance(this.analysisData[h], this.clusters[c]);

                    if (top < 1.0)
                    {
                        top = Epsilon;
                    }

                    // Bottom is the sum of distances from this data point to all clusters.
                    double sumTerms = 0.0;

                    for (int ck = 0; ck < Clusters.Length; ck++)
                    {
                        double thisDistance = CalculateDistance(this.analysisData[h], Clusters[ck]);

                        if (thisDistance < 1.0)
                        {
                            thisDistance = Epsilon;
                        }

                        sumTerms += Math.Pow(top / thisDistance, 2.0 / (this.fuzzynessFactor - 1.0));
                    }

                    membershipMatrix[h, c] = 1.0 / sumTerms;
                }
            }

            this.RecalculateClusterIndexes();
        }

        /// <summary>
        /// Calculates the distance of two vectors.
        /// </summary>
        /// <param name="a">The first vector.</param>
        /// <param name="b">The second vector.</param>
        /// <returns>The calculated distance.</returns>
        private double CalculateDistance(ClusterAnalysisItem a, Cluster b)
        {
            switch (this.distanceMeasureMode)
            {
                case DistanceCalculationMode.Cosine:
                    return DistanceComputation.CosineDistance(a, b.Centroid);

                case DistanceCalculationMode.Covariance:
                    return DistanceComputation.CovarianceDistance(a, b.Centroid);

                case DistanceCalculationMode.DotProduct:
                    return DistanceComputation.DotProductDistance(a, b.Centroid);

                case DistanceCalculationMode.Euclidean:
                    return DistanceComputation.EuclideanDistance(a, b.Centroid);

                case DistanceCalculationMode.Manhattan:
                    return DistanceComputation.ManhattanDistance(a, b.Centroid);
            }

            return 0;
        }

        /// <summary>
        /// Calculate the objective function.
        /// </summary>
        /// <returns>The objective function as double value.</returns>
        private double CalculateObjectiveFunction()
        {
            double jk = 0;

            for (int i = 0; i < this.analysisData.Length; i++)
            {
                for (int j = 0; j < this.Clusters.Length; j++)
                {
                    jk += Math.Pow(membershipMatrix[i, j], this.fuzzynessFactor) *
                          Math.Pow(this.CalculateDistance(this.analysisData[i], this.clusters[j]), 2);
                }
            }

            return jk;
        }

        /// <summary>
        /// Calculates the centroids of the clusters 
        /// </summary>
        private void CalculateClusterCenters()
        {
            for (int j = 0; j < this.Clusters.Length; j++)
            {
                Cluster c = this.Clusters[j];
                double l = 0.0;
                double[] uC = new double[c.Centroid.Count];

                for (int i = 0; i < this.analysisData.Length; i++)
                {
                    double uu = Math.Pow(membershipMatrix[i, j], this.fuzzynessFactor);

                    for (int k = 0; k < c.Centroid.Count; k++)
                    {
                        uC[k] += uu * c.Centroid[k];
                    }

                    l += uu;
                }

                for (int k = 0; k < c.Centroid.Count; k++)
                {
                    c.Centroid[k] = ((int) (uC[k] / l));
                }
            }
        }

        /// <summary>
        /// Executes the cluster analysis.
        /// </summary>
        public override void ExecuteClusterAnalysis()
        {
            this.ExecuteClusterAnalysis(1e-5);
        }

        /// <summary>
        /// Executes the cluster analysis.
        /// </summary>
        /// <param name="precision">The precision of the cluster analysis.</param>
        public override void ExecuteClusterAnalysis(double precision)
        {
            int i = 0;
            const int maxIterations = 20;

            do
            {
                i++;
                this.objectiveFunction = this.CalculateObjectiveFunction();
                this.CalculateClusterCenters();
                this.Step();
                double Jnew = this.CalculateObjectiveFunction();

                if (Math.Abs(this.objectiveFunction - Jnew) < precision)
                {
                    break;
                }
            } while (maxIterations > i);

            this.neededIterations = i;

            for (int j = 0; j < this.analysisData.Length; j++)
            {
                this.clusters[(int) this.analysisData[j].ClusterIndex].ClusterMembers.Add(this.analysisData[j]);
            }
        }
    }
}